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PURPOSE:  This  Coastal  and  Hydraulics  Engineering  Technical  Note  (CHETN)  investigates  the 
ability  of  a  modeling  scheme  to  compute  not  only  the  flow  variables  (pressure  and  velocity),  but 
also  the  location  of  the  free  surface.  The  free-surface  scheme  is  implemented  in  the  three- 
dimensional  (3-D)  Navier-Stokes  module  of  the  Adaptive  Hydraulics  (ADH)  code  and  provides  a 
means  to  model  free-surface  flow  over  a  hydraulic  structure.  A  review  of  existing  methods  of  free- 
surface  tracking  is  provided.  A  particular  form  of  a  moving-mesh  algorithm  is  tested  by  comparing 
simulation  results  of  flow  over  a  weir  with  published  data. 

BACKGROUND:  In  modeling  the  hydrodynamics  of  water  bodies,  the  water  surface  can  be 
regarded  as  a  free  surface  or  as  a  rigid  lid.  Treatment  of  the  free  surface  is  discussed  in  the 
following  paragraphs.  Treating  the  surface  as  a  rigid  lid  (solid  boundary)  requires  that  the  pressure 
be  computed  at  this  boundary.  Modeling  the  water  surface  as  a  rigid  plane  is  reasonable  for  large, 
deep  bodies  of  water  where  the  Froude  number  is  small;  but  it  is  usually  not  appropriate  for  flow 
over  hydraulic  structures  where  significant  variations  in  the  water  surface  are  expected. 

Calculation  of  the  free  surface  in  3-D  Navier-Stokes  models  may  be  complicated,  depending  on 
the  hydraulic  conditions  of  the  problem  at  hand.  Aliabadi  et  al.  (2003)  review  methods  used  to 
address  the  free-surface  issue  in  multidimensional  modeling.  They  report  two  methods  for 
representing  a  moving  free  surface:  the  fixed-mesh  technique  and  the  moving-mesh  technique.  The 
choice  of  technique  depends  on  the  complexity  of  the  expected  free-surface  shape.  The  question  is 
whether  small  water-surface  displacement  is  expected,  or  whether  breaking  waves  and  hydraulic 
jumps  are  expected.  Moving-mesh  algorithms  employ  a  kinematic  boundary  condition  on  the  free 
surface  such  that  the  mesh  moves  with  the  free  surface.  A  commonly  used  technique  is  the 
Arbitrary  Lagrangian-Eulerian  (ALE)  description  (Donea  1982;  Farhat  et  al.  1995).  Another 
moving-mesh  method  uses  space-time  finite  element  formulations  of  the  governing  equations 
(Gulcat  et  al.  1997;  Tezduyar  et  al.  1992).  Behr  (2001)  used  the  space-time  finite  element 
approach  to  capture  the  contorted,  time-dependent  free  surface  of  flow  passing  over  a  fixed-crest 
weir.  Fixed  mesh  methods  include  the  popular  volume  of  fluid  (VoF)  (Hirt  and  Nichols  198 1),  the 
level  set  method  (LSM)  (Osher  and  Sethian  1988;  Sussman  et  al.  1994;  Sethian  1999;  Chang  et  al. 
1996;  and  Lin  et  al.  2005),  and  the  interface-sharpening/global  mass  conservation  (IS-GMC) 
methods  (Aliabadi  et  al.  2003). 

The  VOF  method  and  LSM  do  not  determine  the  water  surface  position  explicitly.  The  primary 
advantage  of  these  fixed  mesh  methods  is  that  the  water  surface  need  not  be  smooth  or  even 
single-valued.  These  methods  are  robust  even  for  their  contorted  surface  patterns  including  wave 
breaking.  Moving  mesh  methods  can  be  posed  to  provide  accurate  water-surface  location,  but  at 
the  risk  of  mesh  tangling  (or  element  shearing).  The  recent  paper  by  Lin  et  al.  (2005)  illustrates  the 
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applicability  of  the  LSM  to  mechanical  and  civil  engineering  applications.  They  show  that  LSM 
can  be  customized  for  solution  of  flow  over  hydraulic  structures. 

Aliabadi  et  al.  (2003)  have  combined  moving-mesh  and  fixed-mesh  techniques  in  a  finite  element 
context  to  simulate  the  physical  motion  of  moving  objects  influenced  by  the  flow  field,  leading  to 
a  coupling  of  fluid/structure  interaction.  For  surface  evolution  they  use  the  interface-sharpening/ 
global  mass  conservation  (IS-GMC)  approach,  which  is  a  variant  of  the  VoF  approach  designed 
primarily  to  avoid  the  somewhat  diffuse  phase  interfaces  that  simple  VoF  methods  tend  to  produce 
while  maintaining  the  discrete  mass  conservation  properties  characteristic  of  VoF. 

APPROACH:  The  3-D  Navier-Stokes  module  of  the  ADH  flow  solver  is  used  to  model  the 
velocities,  pressures,  and  water-surface  elevations  of  the  flow  approaching  and  overtopping  a 
removable  spillway  weir  (RSW).  ADH  is  a  suite  of  finite  element  models  of  which,  the  Navier- 
Stokes  module  solves  the  Reynolds-Averaged  Navier-Stokes  (RANS)  equations  using  the  method 
of  finite  elements.  ADH  has  been  applied  to  flows  approaching  the  Ice  Harbor  and  Lower 
Monumental  Dams  (Stockstill  et  al.  2005a,  2005b). 

The  3-D  numerical  flow  model  reproduced  a  l:25-scale  physical  model  studied  by  ENSR 
International  (2005).  The  physical  model  is  a  section  model  of  the  Lower  Monumental  Dam 
having  three  spillway  bays.  The  RSW  was  placed  in  the  center  spillway  bay.  The  computational 
issue  addressed  in  this  study  is  whether  ADH  can  successfully  simulate  the  free-surface  flow 
approaching  and  overtopping  the  RSW.  This  test  compares  a  solution  from  the  numerical  model  to 
observed  laboratory  data. 

The  ADH  Navier-Stokes  code  employs  a  moving-mesh  method  similar  to  that  described  by  Behr 
(2001).  ADH  uses  the  space-time  finite  element  approach  to  capture  the  time-dependent  free 
surface.  This  is  essentially  an  ALE  method  of  describing  the  fluid  dynamics.  The  mesh  nodes  are 
restricted  to  move  only  in  the  vertical  direction.  Given  a  specified  pressure  at  the  free  surface 
(typically  zero),  the  kinematic  boundary  condition  is  satisfied  and  the  displacement  of  free-surface 
nodes  is  computed.  That  is,  the  unknown  at  free-surface  nodes  is  displacement  rather  than 
pressure.  A  displacement  of  the  free  surface  is  computed  each  time-step.  To  help  preserve  mesh 
quality,  nodes  on  the  interior  of  the  domain  are  moved  (vertically)  as  a  function  of  the  free-surface 
displacement.  Certain  nodes  are  chosen  as  fixed  thus  preserving  the  geometry  of  the  dam  (lower 
portions  of  the  RSW,  the  RSW  crest,  spillway  pier  bases,  etc.). 

Steady  state  is  established  from  an  initial  condition  of  a  quiescent  pool  and  an  assumed  pool 
elevation.  Flow  is  first  established  with  a  fixed  mesh  (rigid  lid)  then  the  mesh  is  allowed  to  move 
and  the  model  is  run  further  in  time  until  steady  state  with  the  free  surface  conforming  to  the 
condition  of  zero  pressure  on  the  water  surface. 

The  flow  boundaries,  including  the  geometric  features  of  the  spillway  piers  and  the  RSW,  were 
constructed  in  a  computer-aided  design  (CAD)  model.  This  CAD  description  of  the  geometry  was 
then  used  as  input  for  the  grid  generator.  The  3-D  volume  mesh  is  of  the  unstructured  type  using 
tetrahedral  elements.  Figures  1  and  2  are  scenes  of  the  surface  mesh  making  up  the  model.  Three 
spillway  bays  were  reproduced  to  simulate  the  physical  model.  The  computational  mesh  developed 
to  model  the  l:25-scale  sectional  model  consisted  of  78,300  nodes  and  412,668  tetrahedral 
elements.  The  numerical  model  was  run  at  laboratory  scale  ( 1 :25),  but  the  results  in  this  report  are 
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given  in  prototype  units  as  they  are  in  the  ENSR  report  (ENSR  International  2005).  The  surface 
meshes  shown  in  Figures  1  and  2  partially  illustrate  the  resolution  of  the  computational  mesh.  The 
surface  mesh  is  composed  of  the  individual  faces  of  the  tetrahedral  elements  that  form  the 
boundary. 

A  forebay  water-surface  elevation  of  164.6  m  (540  ft)  and  a  discharge  of  265.6  cu  m/sec 
(9,375  cfs)  were  simulated.  The  boundary  conditions  consist  of  setting  flow  rates  at  the  inflow,  at 
the  outflow,  at  the  water  surface,  and  at  all  solid  boundaries  comprising  the  computational  domain. 
Boundary  conditions  at  the  solid  boundaries  prohibit  the  flux  of  mass  and  momentum  through  the 
boundary.  This  no-flux  boundary  condition  is  applied  to  the  solid  boundaries  and  a  drag  coefficient 
is  assigned  to  the  faces  of  each  element  forming  the  solid  boundaries.  The  water  surface  is 
modeled  as  a  free  surface  that  moves  in  space.  This  requires  that  the  displacement  be  computed 
and  that  the  pressure  at  nodes  on  the  water  surface  be  specified.  A  stress-free  condition  is  applied 
for  the  tangential  component  of  shear  at  the  water  surface. 

The  average  velocity  was  specified  uniformly  across  the  upstream  boundary.  Nodes  along  the 
water  surface  at  the  inflow  boundary  were  fixed  and  the  pressure  at  these  nodes  was  set  to  zero. 
Outflow  discharge  was  set  as  natural  boundary  condition  through  each  outflow  face.  Nodes  along 
the  water  surface  at  the  outflow  were  allowed  to  move.  The  pressure  at  each  node  on  the  water 
surface  was  specified  as  zero  when  they  were  allowed  to  move. 

RESULTS:  Comparisons  of  the  ADH  results  and  the  laboratory  data  are  provided  in  Figures  3-6. 
The  observed  laboratory  data  are  labeled  as  ENSR  data  in  the  figures  and  computed  results  are 
labeled  ADH.  Water-surface  comparisons  are  presented  in  Figure  3.  The  profile  along  the  center 
line  shows  that  the  computed  drawdown  is  close  to  that  measured  in  the  physical  model.  Figure  4 
is  a  plot  of  the  observed  and  computed  pressure  along  the  face  (bed)  of  the  RSW.  These  pressures 
are  presented  as  pressure  head  elevations.  The  numerical  model  accurately  captures  the  pressure 
variations  as  the  flow  approaches  the  RSW  crest.  This  is  where  the  flow  lines  are  converging 
(Figure  7).  However,  after  the  flow  passes  over  the  crest,  the  streamlines  appear  to  diverge 
(Figure  7),  and  the  computed  pressure  magnitudes  are  greater  than  those  measured.  This  might  be 
attributed  to  the  choice  of  eddy  viscosity  or  the  fact  that  the  outflow  boundary  is  too  close  to  this 
area.  The  outflow  boundary  location  was  chosen  so  that  flow  approaching  and  overpassing  the 
RSW  crest  could  be  accurately  modeled.  The  accurate  modeling  of  pressures  on  the  RSW  face 
downstream  of  the  crest  was  not  considered  as  important  for  this  study  as  the  accurate 
reproduction  of  the  velocity  approaching  the  weir. 

The  computed  and  the  measured  velocities  can  be  compared  by  examining  the  velocity  vectors 
plotted  in  Figure  5.  In  general,  the  comparisons  are  favorable  and  show  the  computational  model 
does  a  good  job  of  computing  the  velocity  vectors.  Observed  and  computed  velocity  magnitude 
contours  shown  in  Figure  6  illustrate  the  velocity  distribution  along  a  vertical  plane  placed 
perpendicular  to  the  axis  of  the  dam  at  the  center  of  the  RSW.  The  solution  is  shown  in  the  scene 
provided  in  Figure  8,  which  is  an  oblique  view  of  the  flow  approaching  and  overtopping  the  RSW. 
The  ceiling  of  the  domain  was  found  by  the  moving  mesh  which  conforms  to  the  shape  of  the  free 
surface.  The  stream  lines  show  that  the  flow  is  3-D,  having  both  lateral  and  vertical  components  as 
a  particle  moves  in  the  longitudinal  direction  over  the  weir. 
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CONCLUSIONS:  This  CHETN  describes  the  findings  of  a  3-D  numerical  model  study  of  flow 
approaching  and  passing  over  a  weir.  The  Navier-Stokes  module  of  the  ADH  code  was  employed 
with  a  moving-mesh  algorithm  to  capture  the  position  of  the  free  surface  as  flow  overtops  the  weir. 
Computational  model  results  are  compared  to  data  obtained  on  a  physical  model.  The  comparisons 
show  that  the  3-D  moving-mesh  model  reproduces  the  flow  approaching  and  overtopping  a  weir 
with  reasonable  accuracy. 

FUTURE  WORK:  The  water-surface  displacement  in  this  example  is  both  smooth  and  single¬ 
valued.  Nevertheless,  several  weaknesses  in  the  simulator  were  discovered  through  numerical 
experimentation,  including  outright  failures  of  the  nonlinear  solver  algorithm  under  certain  choices 
of  time-step,  preconditioner,  and  mesh  parameters.  Jacobian  matrices  in  finite  element 
computations  can  be  ill-conditioned  due  simply  to  poor  mesh  quality,  and  this  ill-conditioning 
leads  in  turn  to  dynamic  instability.  Numerical  models  relying  on  explicit  parameterizations  of  a 
free  surface  (e.g.,  moving-mesh  and  ALE  formulations)  are  prone  to  develop  this  type  of  mesh- 
induced  instability  (Sethian  1999).  For  this  reason,  and  because  even  more  complex  free  surface 
and  fluid-structure  problems  are  widespread,  the  ADH  team  is  developing,  independently  of  this 
investigation,  a  set  of  robust  and  accurate  tools  for  implicit  moving  boundary  problems.  These 
methods,  along  with  hybrid  methods  that  can  maintain  a  high  quality  explicit  surface  mesh,  will 
reduce  the  effort  and  increase  the  accuracy  for  simulations  such  as  this  one. 

POINTS  OF  CONTACT :  This  CHETN  was  developed  within  the  Navigation  Systems  Research 
and  Development  Program  administered  by  James  E.  Clausner,  Coastal  and  Hydraulics 
Laboratory,  U.S.  Army  Engineer  Research  and  Development  Center,  3909  Halls  Ferry  Road, 
Vicksburg,  MS  39180.  For  additional  information  contact  Drs.  Richard  L.  Stockstill 
( Richard. L. Stockstill@erdc.usace.army.mil ),  Christopher  Kees  ( Christopher. E. Kees@erdc. usace. 
army.mil),  or  Charlie  Berger  ( Charlie.R.Berger@erdc.usace.army.mil ),  all  of  the  Coastal  and 
Hydraulics  Laboratory. 

This  CHETN  should  be  cited  as  follows: 

Stockstill,  R.  L.,  C.  E.  Kees,  and  R.  C.  Berger.  2006.  Modeling free- 
surface  flow  over  a  weir.  ERDC/CHL  CHETN-XIII-1,  Vicksburg, 

MS:  U.S.  Army  Engineer  Research  and  Development  Center. 

http://ch  l.  erdc.  usace.  army,  mil/chetn/. 
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Figure  1 .  Overall  view  of  surface  mesh  of  Lower  Monumental  section  model  (inflow  and 
water-surface  boundaries  removed). 


Figure  2.  Close-up  view  of  surface  mesh  of  RSW  (outflow  and  water-surface  boundary 
removed). 
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ENSR  Velocity  Magnitude  Contours 
and  Velocity  Vectors 


Figure  6.  Comparison  of  computed  (ADH)  and  observed  (ENSR)  velocity  magnitude  contours. 
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Figure  8.  Velocity  magnitude  contours  and  stream  traces  of  flow  over  RSW. 
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